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This letter confirms the existence of heat transfer enhancement (HTE) and reduction (HTR) in 
turbulent natural convection with polymer additives. HTE and HTR were numerically predicted by 
Benzi et aZ.(PRL, 104 024502, 2010) in homogenous turbulent convection, but experiments by Ahlers 
& Nikolaenko(PRL, 104 034503, 2010) in turbulent natural convection observed HTR only. Using 
direct numerical simulation of natural convection, the present study reconciles earlier numerical and 
experimental work on the basis of the dominant role of polymer length in the polymer dynamics in 
extensional flows. 

PACS numbers: Valid PACS appear here 



In the turbulence research community, high-molecular 
weight polymer additives are well known for their ability 
to reduce turbulent dragp] [2]. Their effect on turbu- 
lent heat transport is far less documented, even though 
polymer solutions may be of interest in the management 
of heat transfer via turbulence control. A recent exper- 
imental study[3] demonstrated that dilute solutions of 
polymers have the ability to produce heat transfer re- 
duction (HTR). A direct numerical simulation, and shell 
model simulation of homogenous natural convection ^ 
at Pr = 1 reported a non-monotonic behavior of heat 
transfer as a function of the Weissenberg number We, 
ratio of the relaxation time scale of the polymer solu- 
tion to the turbulent time scale of the flow. Using direct 
numerical simulation (DNS), heat transfer enhancement 
(HTE) was observed for all simulated We < 1 calculated 
from the flow time scale based on the rms of velocity fluc- 
tuations. Benzi et al. also used a shell model to expand 
the range of We and predicted HTR at large We. Ahlers 
& Nikolaenko[3] speculated that HTE might be the result 
of the absence of walls in the simulation. Hereafter the 
modification caused by polymer addition to heat trans- 
fer, measured by the Nusselt number Nu the ratio of 
the convective to conductive heat fluxes, is defined as 
HTE or HTR = (Nu p /Nu B - 1) x 100 (%) with respect 
to the Newtonian solvent heat transfer under the same 
temperature conditions 

This letter reconciles Ahlers & Nikolacnko[3 's HTR 
measurements and Benzi et a/.[4j's predictions of HTE 
and HTR, by using direct numerical simulations of nat- 
ural convection in a polymeric fluid between two infinite 
horizontal, isothermal walls. Additionally, the present 
study demonstrates the importance of the polymer length 
L in the selection of the regime of heat transfer (HTE or 
HTR), and identifies the specific polymer/flow interac- 
tions that lead to HTE or HTR. 

Turbulent natural convection flows are simulated in 
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a cartesian domain defined by the orthonormal vector 
base (e x ,e y ,e z ) where x, y and z are the two horizontal 
directions and the vertical direction, respectively. The 
components of the velocity vector u are u, v, and w and 
are normalized by the free-fall velocity, or convection ve- 
locity U c — \/agHA, where a, g, H and A the fluid's 
coefficient expansion, gravity, distance between the two 
no-slip, isothermal, horizontal walls, and the tempera- 
ture difference between the bottom (hot) and top (cold) 
walls. In all simulations, H = 1 and A = 1. The 
Rayleigh number, ratio of buoyancy forces to thermal and 
momentum diffusive forces, Ra = agH 3 A/ \vn) is low, 
Ra = 10 5 , yet within the turbulent regime. The kine- 
matic viscosity and thermal diffusivity of the solvent are 
identified as v and k, respectively. The Prandtl number 
is also fixed, Pr — vj n — 7 which is realistic for water. 
The Reynolds number based on U c and H is therefore 
Re = Ra^^Pr^ 1 / 2 = 119. The flow is incompressible 
(V • u = 0), periodic in horizontal directions and the 
buoyancy effect is simulated using the Boussinesq ap- 
proximation, with the following transport equations for 
velocity, pressure p and temperature 9 and temperature 
fluctuations 9' around the hydrostatic temperature pro- 
file : 



du 



+ (u- V)u = -Vp + 9'e z 
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The parameter /? is the ratio of solvent viscosity to 
the zero-shear viscosity of the polymer solution and af- 
fects both the viscous stress and polymer stress terms in 
Eq. ([!]). The polymer stress tensor T is computed us- 
ing the FENE-P (Finite Elastic Non-linear Extensibility- 
Peterlin) model: 



(3) 
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where the tensor C is the local conformation tensor of 
the polymer solution and I is the unit tensor. The prop- 
erties of the polymer solution are f3, the relaxation time, 
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FIG. 1: Modification of the viscoelastic flow Nusselt number 
normalized by the solvent Nusselt number Nu a . Main graph: 
Nu/Nus as a function of the Weissenberg number We for L = 
10 (v), L = 25 (□), L — 50 (A), and L = 100 (CO- Power 

law fits of Nu/Nu s ( ): L = 50, oc We~° - 1 ; L = 100, 

oc We~ ' 2 . Insert: Nu/Nu a as a function of the polymer 
length L for We = 10 (O) and We = 40 (A). 

here based on the convection scales (We = XU C /H) and 
the maximum polymer extension L. The FENE-P model 
assumes that polymers may be represented by a pair of 
beads connected by a non-linear spring and denned by 
the end-to-end vector q. The conformation tensor is the 
phase-average of the tensorial product of the end-to-end 
vector q with itself, C = (q ® q) whose transport equa- 
tion is 

8C 

— +(u- V)C = C(Vu) + (Vu) T C-T (4) 

On the lhs. of Eq. Q, the first two terms are respon- 
sible for the stretching of polymers by hydrodynamic 
forces, while the third term models the internal energy 
that tends to bring stretch polymers to their least en- 
ergetic state (coiled). The FENE-P model has demon- 
strated its ability to capture the physics of polymer drag 
reduction [21 [5]. In the present work, we make the further 
assumption that the thermal conductivity is independent 
of the polymer concentration, which should be reasonable 
within the dilute approximation. 

Eqs. ([l|4| are solved using finite differences on a stag- 
gered grid, following [5]. The code has been validated 
against existing databases of turbulent channel flows 
and natural convection simulations [5J. All simulations 
are performed in a computational domain of dimensions 
8H x 8H x H and resolution 128 x 128 x 129. The 
Newtonian, HTE= +10% (L = 25. We = 10), and 
HTR= —30% were repeated on domains with twice the 
spatial resolution and also doubled lateral dimensions 
with virtually no change in heat transfer. 

In Fig. [T] several simulations were used to map the 
heat transfer response of a wide range of relaxation time 
0.1 < We < 45, and polymer lengths L = 10,25,50,100 
with a few additional simulations for L = 15,20,30,40 
at We = 10 and L = 25,40 at We = 40. Simulations at 
higher We were found to require smaller time and space 



resolutions and will be discussed in subsequent publica- 
tions. As predicted by the shell model of Benzi et al. [4], 
HTE is observed at low We for nearly all polymer lengths 
(except for L — 100), yet the present simulated enhance- 
ment (HTE= +12% for L = 25 and We = 10) is much 
more modest than the maximum observed in Benzi et 
al.'s shell model of homogenous convection, Nu/Nu s w 6, 
or their DNS where Nu/Nu s w 2.4 with L = 30. The 
absence of walls appears therefore to magnify the HTE 
ability of viscoelastic thermal convection flows. Fig. [I] 
also highlights the critical role of the polymer length on 
the heat transfer performance. Long polymers, L > 50, 
show hardly any HTE. For the range of We considered, 
the short polymer simulations show a decrease of the 
HTE effect, yet do not reach HTR for the range of con- 
sidered We. The two HTR flows exhibit a power-law 
behavior at high We, with Nu oc We -0 - 1 and We" - 2 , 
for L = 50 and 100, respectively. Although this result 
indicates a dependence of the power law exponent in L, 
the limited available data does not allow for a definitive 
conclusion. Combining the shell model and the scaling 
theory of Grossmann and Lohse[l7!, Benzi et aL[l] pre- 
dicted the HTR behavior to Nu oc We" 1 . Since the 
shell model prediction were performed at much higher 
Rayleigh numbers, it may be assumed that Benzi's work 
may describe an asymptotic behavior of HTR. 

The topological differences between Newtonian, HTE 
and HTR flows are depicted in Fig. [2] The shape of 
convection cells is identified by the ridges observed in 
an instantaneous isotherm at 9 = 0.85. For the New- 
tonian flow, these ridges are highlighted by contours of 
A* = max|A[|, A^ is the real part of the eigenvalues of 
the velocity gradient tensor Vu. Terrapon et aL[9] intro- 
duced A* as a measure of the ability of turbulent flows 
to stretch polymer molecules. It should be noted that 
regions of large A* are also regions or uni- or biaxial- 
extensional flows[10]. For viscoelastic flows, the isotherm 
is colored by the local polymer stretch, calculated as 
tr(C)/i 2 . The topology of convection cells is dramat- 
ically modified by the presence of the polymers, with the 
emergence of highly organized convection cells. The typi- 
cal horizontal length scale of convection cells is measured 
using the first negative minimum of the radial correlation 
function of vertical velocity fluctuations in planes parallel 
to the walls (not shown). This length scale drops from 
3.2 in the Newtonian flow to 0.87 in the HTR flow and 
0.7 in the HTE flow. The structure of regions of rotation- 
dominated flows is shown by positive isosurfaces of the 
second invariant Q of the velocity gradient tensor [5). The 
Q-criterion is sensitive to the intensity of vortices and 
is adjusted to identify vortical structures in the core of 
the Newtonian flows with Q = 0.1, arising from shear 
layer instability in plumes. Such small scales vortices 
are absent in the HTE flow, and, with a much lower 
threshold Q = 0.01, the structure of convection cells is 
isolated. The vertical slice of polymer stretch contours 
shown in Fig. [2] indicates sustained stretching over most 
of the plumes' vertical extent. The absence of core tur- 
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FIG. 2: Juxtaposed snapshots of domain halves of the 
Newtonian flow (—8 < y/H < 0) and viscoelastic flows 
(0 < y/H < +8) for HTE= +11% (top) and HTR= -30% 
(bottom), looking down to the hot plate at a 45° angle. The 
geometry of the convection cells is highlighted by an isosurface 
of temperature represents 6 — 0.85 colored by the maximum 
magnitude of eigenvalues of the local velocity gradient tensor 
A* (Newtonian flow) or the normalized trace of the confor- 
mation tensor tr(C)/L 2 (viscoelastic flows). The panels at 
x/H = —8 are vertical contour slices of the flow domain using 
A* and tr(C)/L 2 . White isosurfaces of the second invariant 
of the velocity gradient tensor Q show the topology regions 
dominated by rotation[8] over < x/H < 8: Newtonian, 
Q — 0.1; viscoelastic, Q — 0.01. 



bulence is not surprising since polymers have been shown 
to stabilize shear layers [TT]. The HTR flow does not ex- 
hibit as much coherence as HTE. In fact convection cells 
identified by Q = 0.01 isosurfaces appear to break down 
in smaller structures. The vertical slices show a much 
rapid extinction of polyerm stretch away from the ori- 
gin of the plume, as well as more horizontal motion of 
polymer stretch. Using two-dimensional simulations, the 
present author has isolated in plumes a polymer-driven 
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FIG. 3: Profiles of rms of velocity and temperature fluctua- 
tions in the lower half of the computational domain. Newto- 
nian, HTE and HTR flows are denoted by lines without sym- 
bols, closed and open symbols, respectively. u r m S : , A; 



instability at large We and L, which is the likely cause 
of the observed small scales. Since the initial study sug- 
gests that this instability is not a major component of 
the HTR mechansim, this matter will be discussed in a 
future publication. Lastly, Fig.[2]shows much lower poly- 
mer stretch in the boundary layers of HTE convection 
cells than HTR. 

Fig. [3] plots the vertical distributions of the rms of 
temperature, horizontal and streamwise velocity fluctua- 
tions. The respective thickness of the momentum 5 U and 
thermal 5g boundary layers is estimated by the altitude 
of the maxima of M rm s and 9 rms [5] . Newtonian and HTE 
temperature and vertical velocity fluctuations are very 
close, especially in the boundary layer region. The main 
difference is observed in the horizontal velocity, with a re- 
duction of 5 U (S^ TE a 0.12 < 5 e ~ 0.14 < 5™ cwt w 0.18) 
and a reduction of velocity magnitude. The HTR flow 
significantly departs from the Newtonian flow for all 
quantities shown in Fig. [3] The momentum thickness 
is virtually identical to that of HTE but the thermal 
boundary layer is much thicker (Sg ~ 0.22). Temper- 
ature fluctuations are increased in plumes but velocity 
fluctuations are significantly reduced in both plumes and 
boundary layers. The absence of core turbulence in HTR 
and HTE flows leads to convection cell-driven flows which 
is reflected in the fact that rarms sa 2 x u rms , or the the 
intensity of plume velocity fluctuations are the sum of ve- 
locity fluctuations in the boundary layers of two adjacent 
convection cells. 

Fig. prelates the heat flux (w9) z = nA/H(Nu - 1) to 
the elastic energy (e p ) 2 = (u • ( (1 — fi)/Re V • T)^ , where 
(■)z denotes the averaging over time and homogenous di- 
rection at altitude z. Following Grossmann & Lohse[7], 
the focus of the discussion is the relation: 

^{Nu - l)RaPr- 2 = (e u bl ) v + (e u pl ) v - (e p bl ) v - (e p pl ) v 

(5) 

where e u is the TKE dissipation rate, (-)y the averaging 
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FIG. 4: Vertical profiles of the elastic energy contribution to 

the Nusselt number as described in Eq. d5J. -(£y) ? : , 

HTE; , HTE. -(££,}*: O, HTE; A7HTR. Vertical heat 

flux (iu0)»: , HTE; , Newtonian; , HTR. 

over time and the entire flow volume, and the subscripts 
bl and pi denotes the respective contribution of boundary 
layers and plumes. For clarity, Fig. [4] plots only — (e^)* 
and — {£ p p i)zi calculated from the horizontal and vertical 
contributions of e p . In the HTE flow, the volume aver- 
age of the plume elastic energy, —(£pi)v = 0.86 domi- 
nates that of boundary layers — (£m)v = 0.17. The HTR 
flow shows the inverse behavior with —{£ p pl )v = 0.32 and 
— (£ui)v = 0.46. The explanation for the different dis- 
tribution of elastic energy is first found in plumes. The 
Newtonian flow, via contours of A* , shows that the base 
of plumes hosts the largest extensional velocity gradients, 
and polymeric fluids indeed exhibit the most stretch in 
these regions (Fig. |2]). As polymers travel in the plume, 
the velocity gradient drops rapidly, similar to the relax- 
ation of a stretched polymer arising from the sudden 
cancellation of an extensional flow studied by Doyle et 
a/.[12 . These authors demonstrated that the first nor- 
mal stress decreases at a rate of LyJWe/t. Replacing 
time by travel distance in the plume, this decay rate 
can be observes in the vertical slice of the contours of 



tr(C)/L 2 of Fig. [2] showing much more rapid relaxation 
in HTR than HTE. Short polymers are therefore able to 
sustain larger contributions over longer plume lengths. 
The second difference is the significant stretch in bound- 
ary layers at HTR, while almost negligible at HTE. The 
negative (positive) elastic contribution of (ey)* (— (e^)*) 
for z/H + 1/2 > 0.01 is the result of a large increase of 
extensional viscosity driven by accelerating fluid in this 
region. Below the extensional region, boundary layers are 
shear dominated and thus experience the effect of shear 
thinning as shown by positive (negative) elastic contri- 
bution of (e£j)z (— { e \i)z)- Extensional viscosity for a 
FENE-P fluid increases dramatically with increasing L 
and We [15]. which leads to the reduction of horizontal 
velocity fluctuations in the HTR flow (Fig. [3]) . 

In conclusion, the commonality in the present HTE 
and HTR flows is the absence of core turbulence, which 
implies that both flows are driven by convection cells 
alone. The present HTE simulations and the homogenous 
convection simulations of Benzi et al. [4] demonstrate that 
the interactions of polymers with plumes is inherently of 
HTE nature for short polymer length. Long polymers, 
via extensional viscosity in boundary layers and rapid 
relaxation in plumes, cause HTR as observed in Ahlers 
& Nikolaenko[3J's experiment with high molecular weight 
polymers at much higher Rayleigh numbers. The extrap- 
olation of the proposed models to high Reynolds number 
should be valid on the basis that polymers have been 
shown to reduce turbulence intensity and mixing in free 
shear [T4"] and isotropic 1") turbulence, which should be 
the dominant flow patterns in core turbulence of natural 
convection at high Ra. 
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